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ABSTRACT 


A generalized  critical  point 
is  characterized  by  the  vanishing 
of  certain  linear  relationships. 

In  particular,  the  dynamics  near 
such  a point  are  completely  non- 
linear. this  paper*— weyanalyzeS 

fluctuations *at  such  points  of 
spatially  homogeneous  systems. 
WU'disuuMS  thermodynamic  critical 
points Aas  a special  case;  but  the 
main  emphasis  is  on  stochastic 
kinetic  equations.  We-  show  that  . - 
fluctuations  at  a critical  point 
cannot  be  characterized  by  a 
Gaussian  density,  but  more 
sophisticated  densities  yield 
reasonable  results?; ‘keji  theory 
is  applied  to  the  critical 
harmonic  oscillator. 
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UNIFORM  TREATMENT  OF  FLUCTUATIONS  AT  CRITICAL  POINTS 


SECTION  1.  INTRODUCTION 


A generalized  critical  point  can  be  characterized  by  the 
vanishing  of  certain  linear  relationships.  Such  a point  may  be 
thermodynamic  or  kinetic.  For  example,  at  the  liquid-vapor  tran 
sition  temperature,  it  is  well  known  that 


In  general,  if  4>(x)  is  some  generalized  potential  and  <J>(0)  =0 
then  the  Onsager  theory  of  non-equilibrium  processes  indicates  that 
perturbations  from  x = 0 evolve  according  to 


In  the  vicinity  of  equilibrium,  x 


f I t 

At  a generalized  critical  point,  4>"  » o.  Hence  the  dynamics  in  the 
vicinity  of  a critical  point  are  highly  nonlinear.  Namely,  one  must 
use  extra  terms  in  the  Taylor  expansion  (1.3),  and  obtain  non- 
. linear  dynamics  at  the  critical  point. 

In  this  paper,  we  consider  thermodynamic  and  kinetic  critical 
points.  As  the  term  is  usually  used,  "critical  point"  refers  to  a 
point  in  parameter  space  that  is  characterized  by  i)  a slowing 
down  of  the  dynamics  and  ii)  long  range  spatial  correlations.  In 
this  paper,  we  do  not  consider  the  long  range  order,  but  restrict 
the  problem  to  spatially  homogeneous  systems,  so  that  only  critical 

[slowing  down  will  be  evidenced.  These  equilibria  (thermodynamic 
case)  or  steady  states  (kinetic  case)  are  characterized  by  non- 
linear dynamics.  This  point  is  given  much  emphasis  in  the  paper. 

Due  to  the  nonlinear  dynamics  the  analysis  of  systems  near  critical 
points  is  quite  difficult.  Some  analysis  has  been  done,  by  Kubo 
et.  al. (1)  and  Nitzan  et.  al. (2) . in  the  first  paper,  two  types 
of  critical  steady  states  were  defined  for  one  dimensional  systems. 
We  briefly  review  these.  Let  x=0  be  the  steady  state.  The  steady 
state  is  of  the  marginal  type(l)  if  perturbations  from  x=0  behave  as 

• 

x ^ ±x2.  (1. 5) 
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The  origin  is  a critical  type  steady  stated)  if  perturbations 
from  x=0  behave  as 


X 'V'  ±x 


(1.6) 


Recently,  the  terminology  of  Kubo  et.  al.  has  been  generalized  to 
multidimensional  systems (3).  The  paper  of  Nitzan  et.  al. (2)  is 
contained  as  a special  case  of  this  paper  and  the  accompanying 
one (4) . 

In  this  paper,  we  will  analyze  fluctuations  at  critical 

points.  Often  one  reads  that  fluctuations  "become  unbounded"  or 

"grow  anomalously"  at  critical  points.  These  statements  are  meant 

in  the  following  sense:  if  one  tries  to  describe  fluctuations  at 

a critical  point  by  a Gaussian  approximation,  then  the  second 
2 

moment  <x  > is  infinite.  We  will  show  that  the  Gaussian  descrip- 
tion of  fluctuations  implicitly  assumes  linear  dynamics.  Since 
critical  point  dynamics  are  nonlinear,  one  should  not  expect  the 
Gaussian  approximation  to  be  valid.  Hence,  the  anomaly  is  not  in 
the  physics,  but  in  the  improper  use  of  mathematical  approxima- 
tions. We  show  that  although  the  Gaussian  approximation  is  not 
valid,  more  complicated  densities  are  appropriate.  The  technique 
to  demonstrate  this  will  utilize  formal  asymptotic  methods. 

Our  results  are  analogous  to  problems  in  optics  (at  a caustic) 
and  wave  mechanics  (at  a classical  turning  point).  In  those  cases. 


the  geometrical  optics  and  WKB  solutions  break  down,  yielding 


i 


infinite  amplitudes.  In  reality,  the  intensity  of  light  at  a 
caustic  is  not  infinite,  but  is  large(5).  At  a caustic,  geo- 
metrical optics  must  be  replaced  by  Airy  or  Pearcey  integrals (5,6) . 
Similar  analyses  hold  at  the  classical  turning  point  (e.g. (7) , (8) ) . 

In  §2,  we  give  the  uniform  treatment  of  fluctuations  at  the 
thermodynamic  critical  point  of  a homogeneous  system.  As  an 
example,  we  calculate  volume  fluctuations  of  van  der  Waals  gas  at 
the  critical  point.  Our  theory  is  equivalent  to  the  classical 
theory  of  phase  transitions ( 8) . The  main  focus  of  this  paper,  how- 
ever, is  kinetic  critical  points,  which  have  a much  richer  dynamical 
behavior.  In  §3,  we  introduce  the  stochastic  kinetic  equation  and 
diffusion  approximation.  The  theory  given  here  is  a variation 
of  the  mode-mode  coupling  theory(lO).  We  discuss  a possible  resolu- 
tion of  the  present  controversy  regarding  a "proper"  expansion  of 
the  Master  Equation  to  obtain  a Fokker-Planck  equation ( 11) . A small 
parameter  arises  in  the  derivation  of  the  diffusion  approximation; 
it  characterizes  the  intensity  of  fluctuations.  The  fluctuations 
are  characterized  by  a density  that  satisfies  the  forward  or  Fokker- 
Planck  equation.  In  this  paper,  techniques  for  the  construction  of 
solutions  of  the  forward  equation  are  given  when  the  underlying 
deterministic  dynamics  exhibit  critical  behavoir.  In  §4  we  de- 
rive solutions  of  the  one  dimensional  time  invariant  Fokker-Planck 
equation.  We  obtain  an  exact  result,  which  is  then  analyzed 
by  asymptotic  methods.  We  obtain  a Gaussian  density  at  a 
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non-critical  steady  state,  an  Airy  density  at  a marginal  type 
steady  state,  and  a Pearcey  density  at  a critical  type  steady 
L_ate.  In  §5,  these  densities  are  used  in  a general  ansatz  ("ray 
method" (12))  to  provide  asymptotic  solutions  of  the  time  dependent 
multidimensional  Fokker-Planck  equation.  We  construct  densities  in 
which  susceptibilities  (i.e.  first  derivatives)  at  the  critical 
point  are  large,  but  finite.  The  same  result  applies  to  variances. 
In  §6,  we  show  how  our  results  can  be  used  to  construct  time- 
dependent  correlation  functions.  In  §7,  we  discuss  an  example  of 
the  critical  harmonic  oscillator (13)  and  show  how  the  correlation 
function  is  constructed. 

Antecedents  to  this  work  are  found  in  Kubo  et.  al. (1) , 
Kitahara ( 14) , Keizer(15)  and  Nitzan  et.  al. (2) . The  present  work 
generalizes  the  results  of  the  above  papers. 


SECTION  2.  UNIFORM  THERMODYNAMIC  THEORY 


In  this  section,  we  derive  the  uniform  theory  for  thermo- 
dynamic critical  points  of  spatially  homogeneous  systems.  The 
ideas  which  arise  here  are  very  similar  to  the  more  complex  ones 
that  arise  in  the  kinetic  case. 


2.1.  GENERAL  THEORY 

The  thermodynamic  theory  proceeds  from  the  Einstein 
fluctuation-entropy  formula.  We  assume  that  the  entropy  of  the 
system  can  be  characterized  by  a parameter  x.  The  equilibrium 
entropy  is  Sn  = S(x_),  Let 


v(x)dx  = Pr  | system  reaches  a state  in  which 
x€(x,x+dx)|  . 


(2.1) 


Then (9) 


s(x)  ” so  TasI 

v(x)dx  = c exp  ^ dx=c  exp  dx. 


(2.2) 


where  c is  a normalizing  constant.  As  is  usually  done,  we  have 
eliminated  the  time  variable  from  the  discussion  of  thermodynamic 
problems (16 ) . This  elimination  has  certain  conceptual  drawbacks 
when  one  tries  to  describe  the  time  evolution  of  a system.  How- 
ever, we  shall  follow  standard  notation  here.  We  also  follow  the 
standard  procedure  of  dividing  the  universe  into  the  sub-system 
of  interest,  characterized  by  a variable  y and  certain  intensive 


parameters  (ct^,  ...  a^)  = a , and  an  external  reservoir  which  is 
assumed  to  remain  in  thermodynamic  equilibrium (17,  page  274).  Then 

v (y ) dy  = c exp  j dy,  (2.3) 

where  AW(y)  is  the  work  done  on  the  subsystem  by  an  external  source 
Let  $(y,a)  denote  the  potential  of  the  system,  so  that  (with  y=yQ 
denoting  equilibrium) : 

AW (y ) = $(y,a)  - *(yQ,a).  (2.4) 


Then 

v (y ) dy  = c exp  |~~<I>^Y>ot)  j dy.  (2.5) 

Usually  (9,16,17)  4>(y,a)  is  expanded  in  a Taylor  series,  keeping 

terms  of  second  order 

<J>(y,a)  = $(y0,a)  + 4>  • (yQ , a)  (y-yQ) 

(yn»c)  2 

+ <t>"  — ^ (y-yQ)  + ...  (2.6) 

At  equilibrium  4»'(yQ,a)  =;  0 and  <t>"(yg,a)  > 0.  Thus  we  obtain 

v(y)dy  ^ exp  T 4>"  (yn , a)  (y-y^ ) 21  dy. 


(2.7) 


Equation  (2-7)  gives  a locally  Gaussian  den.'ity  with  variance 


2 

a = 


kT 


4>"  (yQ/Ct) 


(2.8) 


The  Onsager-Machlup  theory  of  irreversible  processes  (18)  proceeds 
from  this  point. 

Suppose,  however,  that  there  is  a value  a=am  such  that 


$,(yO,0tm)  = *"  (y"'aJ  = 0 n'aJ  * °* 


O'  m 


0'  m' 


(2.9) 


Such  a point  might  correspond  to  a second  order  phase  transi- 
tion (2,9).  A third  term  is  needed  in  the  Taylor  expansion  (2.6) 

For  a near  a , instead  of  (2.7),  one  obtains 
m 


v(y)dy  ^ exp 


=1  j ,y-y  ,2 

kT  | 2 ^Y  Y0 


*’"(y0,a)  3 

+ s-2— 


dy. 


A simple  change  of  variables  converts  (2.10)  to 


v(y)dy  ^ exp 


(2.10) 


[ {z(y)3  “ «z(y)  4 dz(y)  , 


(2.11) 


where  z(y)  is  a regular  function  of  y and  a (a)  is  a regular  func- 
tion of  a with  the  property  that  a (ot  ) = 0.  A density  of  the  form 

m 

(2.11)  is  called  an  Airy  density. 


$ 


At  thermodynamic  critical  points,  if  <J>"(y0/a)  vanishes,  then 
usually  ' (y0,ci)  also  will  vanish.  This  can  be  shown  by  using 
free  energy  arguments (17) . Such  will  not  be  the  case  for  kinetic 
equations,  however. 

It  is  also  possible  that  at  a different  value  of  a,a=a  . 

c 

4>"(y0,ac)  = $',,(y0,ac)  * 0;  $(lV)(y0,ac)  * °*  (2.12) 


Such  a point  corresponds  to  a first  order  phase  transition (2) . We 

take  another  term  in  the  Taylor  expansion  (2.6)  and  obtain,  for 

a near  a : 

c 


v(y)dy  ^ exp 


zi 

kT 


*lv<y0,o)  (y~yo)4  <t>"'  (yo,a)  (y_yo)3 

24 + 6 


+ 


^"(yo.a)  (y-y0)2 

2 


(2.13) 


Equation  (2.13)  can  be  put  into  the  form 

I zi£ii “ila)z<y)2 

kT 


v(y)dy  ^ exp 


■1  I z (y) 4 , “l<o)z(y>2 

;T  j 4 2 a2 (a)  j 


dz  (y) . 


(2.14) 


In  (2.14),  z(y)  is  a regular  function  of  y;  a^,a2  are  re9ular  func- 
tions of  a and  vanish  at  a=a  . We  call  the  density  (2.14)  a 

c 

Pearcey  density. 
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Equations  (2.10,14)  represent  a formal  extension  of 

Onsager's  theory  to  critical  point  phenomena.  In  light  of  (1.4), 

we  are  lead  to  totally  nonlinear  dynamics  when  a=a  or  a=a  . Our 

c m 

result  is,  of  course,  purely  *formal  and  is  applicable  to  small 
deviations  from  equilibrium  only,  which  is  the  best  that  one  can 
expect  of  a thermodynamic  theory. 


2.2.  VOLUME  FLUCTUATIONS  OF  A VAN  DER  WAALS  GAS 


As  an  example  of  the  above  analysis,  we  consider  the  volume 
fluctuations  of  a gas  at  the  liquid-vapor  critical  point,  using  a 
van  der  Waals  model.  Levich(17)  shows  that  in  this  case 


AW  = PqAV  + AF, 


(2.15) 


where  P_  is  the  equilibrium  pressure  of  the  reservoir.  Expanding 


AF  gives 


AW  = PqAV  + 


(II 


AV  + 


92f\  (AV)2 


, { 33f\  (AV) 3 . ( 92f\  (AV)4 

*WL—  + ' 


(2.16) 


where  AV  = V-V^  is  the  deviation  from  the  equilibrium  volume. 
Usually,  only  two  terms  are  of  (2.16)  are  used.  If  we  set 


• 


P = PQ,  we  find 


10 


(2.17) 


E(<4V»2)  - nwwi?  • 

2 

The  value  of  E (AV)  becomes  infinite  at  T=T  . The  cause  of  the 

c 

divergence  is  purely  mathematical:  only  two  terms  of  the  expansion 
were  used.  At  the  critical  temperature (17) 


Then  (2.16)  becomes 


(2.19) 


Equation  (2.19)  is  general.  We  now  specialize  to  a van  der  Waals 
gas,  for  which 


P = 


RT 

V-b 


(2.20) 


The  conditions  (2.18)  lead  to  the  following  values  of  the  critical 
parameters : 


8a 

">  -7  K D 


T 


V 


3b  P 


a 


(2.21) 


Then 


93p\  _ _v  .0123457a 

,3V3 /T  b5 

c 


Then  we  obtain 


v(AV)  d (AV)  = c exp  ~^> 


d(Av) , 


where  c is  the  normalization  constant.  Hence  we  obtain 


E ((AV)2)  = 


. L *2°*p 


dx 


r 


■m 


r r 

i -yx 


,1/2  f yV^2<  dy 


f ^y4/24  dy 

_ m 


(2.22) 


(2.23) 


(2.24) 


Thus,  the  van  der  Waals  theory  (i.e.  the  classical  theory,  see 
below)  predicts  e((AV)2)  « t3/2.  in  Table  1,  we  give  the  results 
of  calculations  of  E ((AV)  ) for  a number  of  gases.  The  above 
theory  yields  values  of  e((AV)2)  which  are  large,  but  not  infinite. 
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The  thermodynamic  theory  presented  in  this  section  is  a 
"classical"  theory (19)  and  thus  will  not  predict  the  "correct" 
critical  exponents.  This  is  a fault  of  the  use  of  thermodynamic 
theory  per  se.  It  is  not  clear  how  the  thermodynamic  theory  given 
above  could  be  modified  to  yield  the  correct  exponents. In  a recent 
work  Mou  et.  al. (34) , using  the  master  equation,  also  derived 
the  classical  result. 


TABLE  1 

VOLUME  FLUCTUATIONS  AT  THE  CRITICAL  POINT 


V (cc/mole)* 
c 


He 

57.6 

4.416x10 

H2 

65.0 

4.336x10 

N2 

90.0 

4.832x10 

CO 

90.0 

5.060x10 

°2 

74.4 

4.702x10 

C2H4 

127.5 

5.168x10 

• 

C°2 

95.7 

5.147x10 

nh3 

72.4 

6.780x10 

• 

h2° 

45 

1.187x10 

*from  (20,  page  18) 
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((Avr)/vc 

416xl02 

336xl02 

832xl02 

060xl02 

702xl02 

168xl02 

147xl02 

780xl02 

2 
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SECTION  3.  STOCHASTIC  KINETIC  EQUATIONS,  DIFFUSION 
APPROXIMATION  AND  FOKKER-PLANCK  EQUATION 

The  thermodynamic  theory  of  §2  can  not  be  used  to  treat 

highly  nonequilibrium  kinetic  phenomena,  which  are  of  interest  in 

. . % 

many  areas  of  chemistry,  physics  and  biology.  Let  x(t)  denote  the 
statistical  variables.  Often  we  can  postulate  an  equation  for  the 

f\j 

mean  value  of  x(t),  x(t): 


dx  , i , . 

dt"  = b (x,a) 


x1 (0)  = Xq 


i ~ 1 f • • • / n 


a = {a1,  ...  %}• 


(3.1) 


In  order  to  treat  fluctuations,  we  need  to  know  the  kinetic  equa- 
tion  that  x(t)  satisfies.  Ideally,  we  would  start  with  the 
Liouville  equation  and  derive  the  kinetic  equation.  Such  a deriva- 
tion is  possible  for  only  the  simplest  system (33) . Instead,  we 
shall  use  a generalization  of  the  Langevin  method.  We  will  add 
a zero-mean  stochastic  term  to  (3.1).  The  stochastic  function 
$(t)  is  characterized  by  a microscopic  time  scale,  x,  small  com- 


pared to  the  macroscopic  scale  on  which  measurements  are  made. 


t 


where  n is  a small  parameter.  We  will  not  assume  that  Y has  a 
6-correlation  function  and  let 


Ykl  = J E ^k  (s)  (o)^  ds. 


'Xi  'Xj 

We  assume  that  x(t)  = x^tt)  satifies  the  stochastic  kinetic 
equation 


(3.3) 


In  equation  (3.3),  e is  a small  parameter  characterizing  the  size 
of  the  system  and  related  to  the  intensity  of  the  fluctuations 
(1,2,13,15).  Hence  e-+0  corresponds  to  the  thermodynamic  limit. 

The  field  f^j(x)  is  a given  deterministic  field.  Ideally,  one 
would  like  to  calculate  f^  from  basic  principles.  Since  (3.3)  is 
somewhat  ad-hoc,  a prescription  is  needed  for  the  calculation  of 
fj.  (One  such  prescription  is  the  fluctuation-description  theorem. 
Another  is  given  by  Keizer(22)).  As  n^O , x^(t)  ■+  x(t),  a diffusion 
process  (22).  We  set 


u(x)  = E. 


x (0)  = x 


(3.4) 


Then  u(x)  satisfies 


i 9u 


9x19x:) 


i 

+ c e 


(3.5) 


In  (3.5),  we  use  the  convention  that  repeated  indices  are  summed 
from  1 to  n and 


13  _ ^i^i  kl  . Ik, 


a J ! ( Y + Y ) 


(3.6) 


1 

C = Y f 


kl„j  8fl  1 9 _i j 

* r = -7  r a 


k ax1’  4 ax3 


(3.7) 


On  the  other  hand,  if  \ has  a {-correlation  function  (white  noise)  , 
then  we  obtain 


3u  = ea13  a2u  .i  au 

at  2 1 ' 3 771 


(3.8) 


ax  ax- 


ax" 


Numerical  work (3)  indicates  that  if  the  boundaries  are  non-singular, 
then  equation  (3.5)  and  (3.8)  yield  equivalent  solutions  for  e<.l. 

We  also  note  that  equation  (3.3)  is  a stochastic  equation 
with  correlations  and  hence  is  a more  reasonable  representation 
then  a white  noise  equation.  Furthermore,  equation  (3.5)  [or  (3.8)] 
is  derived  rigorously — no  expansion  procedure  is  needed  (compare (11) ) 
Equation  (3.8)  is  the  backward  equation.  Usually,  in  the  physical 
literature  the  forward  or  Fokker-Planck  equation  is  used.  This 
equation  can  not  be  derived  rigorously  •.  No  expansion  procedure 
will  rigorously  give  the  Fokker-Planck  equation ( 1 , 11) . Instead, 
we  will  obtain  the  Fokker-Planck  equation  by  using  the  theory  of 
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SECTION  4.  CANONICAL  DENSITIES 


In  this  section,  we  consider  the  time  independent,  one 
dimensional  Fokker-Planck  equation 


e 


(av) 


xx 


2 


(bv)^  = 0 


(4.1) 


subject  to 

/CO 

v(s)ds  = 1 lim  v(s)  = 0.  (4.2) 

,oo  I s | ->“> 

The  solution  of  (4.1,2)  gives  the  steady  state  (but  not  necessarily 
equilibrium)  density  for  a process  satisfying  (3.3).  Our  results 
are  valid  for  £-*-0  ("thermodynamic  limit")  and  generalize  the 
thermodynamic  results  of  §2.  In  later  sections,  we  generalize  the 
solutions  obtained  here  to  solve  time  dependent,  multidimensional 
problems . 

When  (4.1)  is  integrated  twice  and  (4.2)  is  applied,  we  find 

(4.3) 


(4.4) 


v(x)  = k 


exp  | / H ds)_  • 


where  k is  the  normalization  constant 


= l exp\r  § 


ds  } dx. 
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The  main  contribution  to  (4.3)  comes  from  the  maximum  of  the 


function 


-/> 


(4.5) 


We  now  assume  that  there  is  a steady  state  (i.e.  b(x)=0  has  a 
solution)  , Xq . The  steady  state  is  classified  according  to  its 
dynamic  behavior. 

The  normal  type  steady  state  x is  characterized  by 


b(xQ)  =0,  b ' (xQ)  jL  0 


(4.6) 


We  are  interested  in  stable  steady  states,  so  that  we  assume 
b' (xQ)  < 0.  Thus,  perturbations  from  xQ  decay  exponentially. 
When  <t(x)  is  expanded  about  xQ , we  obtain: 


j -|b'(x  ) | (x-x  ) 
v(x)  ^ k exp  { — 


exp  | 


- |b  ' ( xq ) | (x-xQ) 


(4.7) 


Thus,  we  obtain  a locally  Gaussian  density,  for  small  e.  This 
result  has  also  been  derived  by  Kubo  et.  al.(l)  and  Keizer(23)  by 
feren*-  arguments  . it  is  the.  standard  result  in  the  theory  of 
nonequilibrium  thermodynamics  (18) . 


i 


i 


In  the  marginal  case,  b depends  on  one  parameter  a such 
that  when  a = ac  the  marginal  type  steady  state  satisfies 

b (xo ' ac)  = 0 b 1 (xq ' ac^  = 0 b" (Xq , ac)  / 0.  (4.8) 

The  canonical  dynamics  corresponding  to  the  marginal  case  are  (4) 

' 2 

x = x - a.  (4.9) 


The  flow  of  such  dynamics  is  sketched  in  Fig.  1.  We  need  to  re- 
place the  conditions  (4.2)  by: 


lim  v(s) 

S->  — oo 


i 


v(s)ds  = 1 


(4.10) 


where  x < °°  is  an  end  value  for  x (see  also  3)  . 

Since  b'(x0»ac)  = 0,  the  expansion  used  to  obtain  (4.7)  breaks 
down.  Hence  the  Gaussian  density  breaks  down.  In  particular, 
from  (4.7)  we  have 

E(x2)  « l/|b*  (xQ,a|  . (4.11) 


Thus,  when  a = ac>  the  Gaussian  density  yields  an  infinite 
variance.  This  is,  of  course,  a purely  mathematical  divergence 
and  has  nothing  to  do  with  the  physical  problem.  We  take  an  extra 
term  in  the  expansion  of  4>(x,a)  to  obtain 
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v(x)  ^ k exp 


b ' (Xq , ot)  (x-Xq ) 2 b" (xQ,a) (x-xQ) 3 

eaTxZ)  + Jei 


(4.12) 


A change  of  variables  converts  (4.12)  to  the  Airy  density  (2.11) 

The  critical  type  steady  state  is  characterized  by  two 

parameters,  a,B  such  that  when  a = a and  B = B 

c c 

b(x0*ac'ec)  = b'(x0'acec)  = h’^VV^c1  = 0 


b",(x0'ac'ec)  * °' 


(4.13) 


The  canonical  dynamics  of  a critical  type  dynamical  system  are  (4) 

x = ±x3  + ax  + B.  (4. 14) 

In  this  case  a = B =0.  The  flow  of  (4.14)  is  skeched  in 
c c 

Fig.  1. 

It  is  clear  that  the  Airy  and  Gaussian  densities  both  break 

down  when  a = a and  B = B • In  this  case,  we  take  one  more  term 
c c 

in  the  Taylor  expansion  of  4>(x)  and  obtain  the  Pearcy  density  (2.14). 

The  results  of  this  section  can  be  obtained  by  direct  use  of 
Levinson's  theorem ( 26 ) . It  is  clear  that  :he  Gaussion  approximation 
will  be  valid  whenever  |b ' (Xg , a) /b" (x^ , a) | >>1. 


© 


FIG.  1: 


->«  < 


c*  >0  ) ^ >0 


— •- 


(9 


-y>-  4~  ^ 
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A)  THE  MARGINAL  TYPE  DYNAMICAL  SYSTEM  HAS  TWO  STEADY 
STATES  WHEN  A PARAMETER  a > 0,  ONE  DEGENERATE  STEADY 
STATE  WHEN  a = 0 AND  NO  STEADY  STATES  WHEN  a < 0. 

B)  THE  CRITICAL  TYPE  DYNAMICAL  SYSTEM  HAS  THREE 
STEADY  STATES  WHEN  a > 0,  $ > 0,  ONE  DEGENERATE 
STEADY  STATE  WHEN  a = g = 0 AND  ONE  STEADY  STATE 
WHEN  a > 0,  0 > 0. 
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SECTION  5. 


TIME  DEPENDENT,  MULTIDIMENSIONAL  FOKKER-PLANCK  EQUATION 

In  this  section,  we  construct  regular  solutions  of  the  time 
dependent  Fokker-Planck  equation  (3.5).  Since  Ludwig(5)  has  given 
the  construction  for  the  normal  case,  we  only  consider  marginal 
and  critical  type  steady  states.  Our  goal  is  to  construct  densities 
that  have  finite  second  moments.  Exact  definitions  of  marginal 
and  critical  type  steady  states  in  multidimensional  systems  are 
given  in  Appendix  A. 

5.1.  MARGINAL  TYPE  STEADY  STATE 


We  seek  a solution  of  (3.5)  of  the  form 


v(x,t)  = exp 


=1 

e V 3 


— a 4*  (x,  t)/  £ enzn(x,t) 

-1  n=C 


(5.1) 

.1 


The  form  of  (5.1)  is  a "ray  ansatz"(27).  in  it,  4'(x,t),  a,  and 
the  functions  z° (x, t) , z1 (x, t) . . . are  to  be  determined.  In  practice 
we  are  often  interested  in  just  the  first  term  of  (5.1)»  After 
derivatives  are  evaluated,  terms  are  collected  according  to  powers 
of  e.  We  obtain 


I 


0 = exp  — 


+ (a  - 


[-1  <r-  - /e  {’t + bi,,i 

f2>  S2.  f.y.j  2°  (S  - »2)  + exp  U4  - X»>]  |bjj2° 


+ b1z?  + — ' z°  + (a  - 4'2)4'.Z° 


(5.2) 


z.  + a13  ('-2'F'P.  ¥.z°  + (ot  - V2)  V . . z° 
t ' 1 D ID 


+ 2 (a  - 4'2)4'jz°)  - c 


1 (a  - 4'2)4'iz0)| 


0 exp  uV)  ^ . 


The  leading  coefficient  of  e vanishes  if 


4\  + b14' . + (ct  Y — a1\.'V.  = 0 
t l 2 i d 


(5.3) 


Equation  (5.3)  is  a generalized  eikonal  equation(3).  In  the  rest 
of  this  paper,  we  shall  assume  that  the  initial  data  for  v are 
concentrated  at  a point 


v(x,0)  = 6 (x  - Xq ) 


(5.4) 


At  the  deterministic  steady  states,  N,Q  (see  Fig.  2)  we  expect 
that  dY/dt  = H*t  + b^'l^  =0.  We  set  4*2  = a at  those  points.  The 


node  N should  correspond  to  a local  maximum  for  v(x,t).  Hence 

we  set  4'(N)  = Similar  reasoning  leads  to  4*  (Q)  = - 

The  value  of  a is  still  under termined.  It  can  be  obtained 

by  the  following  iterative  procedure (3) . (If  higher  order  terms 

are  to  be  considered,  then  it  is  necessary  to  expand 

a = ^ e^a,  . In  that  case,  all  of  the  parameters  are  determined 
k=0  K 

in  a manner  analogous  to  the  determination  of  a) . We  start  at 

the  node  N,  where  ¥ = ^ a®  , the  first  estimate  for  a.  Equation  (5.3) 

can  be  solved  by  the  method  of  characteristics (12) .The  characteristic 
equations  are: 

ST  = 1 (5.5a) 


dx  , l , .'v  m2x  ij 
55 - = b + (a  - V )p^a 


df  dx  1,^  Uf2,  ij 

di  = pk  ai"  = -(a  " ^ )a 


PiPj 


ds^  = " ^2Pk  (“'1'al^Pipj)  + b'kpi  + ,("-2  ">y~)  a'kpipj) 


(5.5b) 


(5.5c) 


( 5 . 5d) 


Initial  data  is  given  on  an  ellipsoid  surrounding  N.  As  x-K3,  the 


value  of  a should  approach  - ya  . if  it  does  not,  then  a new 

• 'V/  ( 1) 

estimate  a is  needed.  Elsewhere,  we  have  shown  that  iterates 

% 

of  a can  be  determined  by  using  the  method  of  false  position  and 
that  a is  a regular  function  of  the  deterministic  parameter  a (3) . 


i 


When  N and  Q coalesce  (Fig.  2b) , a = 0.  After  the  annihilation 

f\j 

of  N,Q  (Fig.  2c),  a < 0.  The  stochastic  problem  for  a dynamical 

system  similar  to  Fig.  2c  is  an  old  one,  solved  by  Levinson (28) 

and  Ventcel  and  Friedlin (29) . Consequently,  we  restrict  ourselves 

to  the  dynamic  cases  represented  by  Figs.  2a, b. 

2 

From  (5.5b),  we  see  that  if  'F  = a on  a trajectory,  then 
dx1/ds  = b1,  so  that  the  trajectory  is  a deterministic  trajectory. 
In  this  way,  we  will  be  able  to  estimate  deviations  from  a given 
deterministic  trajectory.  When  (5.3)  is  differentiated  with 
respect  to  x and  evaluated  on  a trajectory  we  find: 


d 4/.  i—  . . 

air  + bVk  * V^iVk  - 0 


k = 1,  2 , . . . n 


(5.6) 


In  (5.6)  the  (-)  sign  corresponds  to  trajectories  that  enter  N, 
the  (+)  sign  to  trajectories  that  enter  Q.  At  either  of  the  steady 
states,  we  obtain 


bVk  - 0 


(5.7) 


k = 1 , . . n 


Equation  (5.7)  can  be  solved  to  yield  values  of  4*  at  N or  Q. 

When  N and  Q coalesce,  so  that  a = 0 and  conditions  (A)  hold, 
it  is  possible  to  show  that  the  4^  can  be  calculated.  For  example, 


consider  the  case  of  only  one  spatial  dimension.  Then  (5.7) 
becomes : 


I 


i.e. 


z°(s)  = z°(0)exp|-  | f(s')ds 


* -JS  f (s' ) ds 


(5.11a) 


When  the  initial  data  is  concentrated  at  a point,  Ludwig  has  shown 
that  the  appropriate  initial  data  for  z is  z°(0)  = constant. 

5.2.  CRITICAL  TYPE  STEADY  STATE 

For  the  critical  type  steady  state,  instead  of  (5.1),  we 
seek  a solution  of  (3.5)  of  the  form 


[_l  (t  4,4  ' nr  ■ &*)  £ enzn(x,t). 


v (x,  t)  = exp  ~ - ev  j|  e“z“(x,t).  (5.12) 

In  this  case,  it  is  possible  to  impose  the  conditions  on  v(x,t)  that 


v-*0  as  | x I -*■» 


, f v(x, 


t)  dx  = 1 


(5.13) 


and  take  all  of  Rn  as  the  domain  of  interest.  Instead  of  (5.3), 
we  obtain 

ij 


^t  + bV  + (V2  - at  - 3)  4'i4'j  = 0. 


(5.14) 


The  value  of  f at  the  deterministic  nodes  and  saddle  Q 

(Fig.  3)  is  determined  in  a manner  analogous  to  the  one  used  in 
§5.1.  The  values  of  the  parameters  are  also  determined  in  a 
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similar  fashion.  It  is  possible  to  show  that  all  constructions 
remain  regular  as  the  steady  states  coalesce (3).  The  function 
z°(x,t)  can  also  be  determined  in  manner  analogous  to  the  previous 
case. 


SECTION  6.  CORRELATION  FUNCTIONS  FOR  CRITICAL  TYPE  SYSTEMS 


In  many  physical  problems,  the  object  of  interest  is  the 
correlation  function 

R(t)  = E |x(t)  x(t  + t)|  (6.1) 

= //xx1  Pr  | x ( t)  € (x , x + dx)  , x(t  + t)£(x\x^  + dx^jdxdx1. 
Since  the  process  in  our  problem  is  assumed  to  be  stationary, 

( a,  i\j  \ 

R(t)  = E <x(0)x(t)>  . Now  we  consider  the  conditional  correlation 

function: 

R (t)  = E | x(0)x(t)  |x(0)  = x ) . (6.2) 

xo  ' u ’ 

If  Vg(Xg)dXg  the  initial  density  for  Xq , then  we  clearly  have 

R(x)  = /Rx  (t) vq (xQ) dxQ.  (6.3) 

However, 

R (x)  = / x'^Pr  I x(t)  ^(x1^1  + dx1)  |x(0)  = xn  l dx1 
Xq  ' u ' 

= / x1  v (x1,t)dx1,  (6.4) 

X0 


where  v (x\x)dx^  was  calculated  in  the  previous  section.  Thus 
X0 


>x0<T>  - ' xl  [-?  (nxl,Tl4)-  ¥■  - Hz(x't) 


dx 

x(0)  = x. 


(6. 


Namely,  we  start  the  ray  calculation  at  x = x^  and  integrate  the 
ray  equations 


= 1 = b1+4'3aljp, 

as 


ds 


(6 . b) 


dV 

ds 


pk  35“  di^  " -(3“1:lPipjPk  * bVi 


3 *'k 

+ ~r  PiPj>' 


(6.7) 


until  s = x.  Thus  (6.5)  can  be  evaluated.  The  full  correlation 
function,  obtained  from  (6.4)  is 


R ( x)  = ff  x exp 


l ( y (x}x) 
~e  V 4 


^ 2 
aV 


z°(!t'T|l  W^V 

x(0)  = xQ 


(6.8) 


In  the  next  section,  we  give  an  example  of  such  a calculation. 
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SECTION  7.  CRITICAL  HARMONIC  OSCILLATOR 


In  this  section,  we  consider  a modified  Duffing  oscillator C30) 
in  contact  with  a heat  bath(13,31).  We  shall  use  a stochastic 
equation  of  the  standard  "mode-mode"  form,  but  initially  will 
indicate  how  a more  general  analysis  would  proceed. 

The  Hamiltonian  of  the  system  is 


H 


k (n) x2 
2 


+ 


4 2 

ax  p_ 

4 2m 


+ $.  (r  ,x)  + 
int  'v* 


(7.1) 


In  (7.1),  <J>.  . (r  ,x)  is  the  interaction  potential  of  the  oscillator 
int  'v* 

with  coordinates  (x,p)  and  heat  bath  with  coordinates  (rx,p^).  The 
last  term  represents  the  kinetic  energy  of  the  heat  bath.  The 
motion  at  the  full  system  (oscillator  + heat  bath)  is  generated  by 


l 

r = 


9H_ 

3p. 


Pi  = - 


3H 

3r3 


i = 1, 


N 


(7.2) 


The  motion  of  the  entire  system  occurs  on  a manifold  in  the  phase 
space  given  by  H = E,  where  E is  the  initial  (i.e.  constant)  energy 
of  the  system.  This  manifold,  M,  will  be  bounded  and  compact. 

We  expect  that  the  full  system  is  ergodic(32).  We  are 

interested  in  a submanifold  of  M,  m\  which  is  the  manifold  of 


34 


(x,p)  coordinates.  A possible  projection  operator  from  M to  M 
is 


■■■ r"' Pr  ■■■  p"Hx'p'° °! 


(7.3) 


Namely,  we  "project"  from  On  M1,  we  assume  that  the  follow- 


ing measure  exists 

P ( t , x , p , A)  = Pr 


j(x(t)  , p ( t JCA  | x 


(0)£ (x,x 


+ dx)  , p(0)€(p,p  + dp) 


)• 


(7.4) 


<\j , 


We  have  introduced  x(t),  p(t)  as  random  variables.  This  is  a 
result  of  the  elimination  of  (^  ,jd  ) from  consideration.  By 

averaging,  we  are  treating  the  latter  variables  as  random;  thus 
x,p  become  random  variables 


Next, 

we  assui 

lim 

At->0 

1 , 

f 

At  J 

1 (x 

f % 

lim 

At->-0 

At  J 

1 (P 

lim 

At+0 

1 i 

At  J 

' (P 

rO, 


/\J  , 


x)  P (t  + At , x, p, dx, dp)  = p/m 

(p  - p)  P (t  + At,x,p,dx,dp)  = 


(7.5) 


-k(n)x  - axJ  - y (x)  p 


.dp)  = Ea, 


(7.6) 

(7.7) 
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and  that  all  other  moments  are  zero.  These  assumptions  have  yet 
to  be  verified  for  any  but  the  simplest  system(33).  These  assump- 


tions lead  to  the  Langevin  equations 


= -k (n) x - 


3 / — dy 

- yx  + yea 


(7.8) 


(7.9) 


where  y(t)  is  an  approximation  to  Gaussian  white  noise  (e.g.(22), 
page  651) . 

The  Fokker-Planck  equation  for  the  density  f(t,x,p)  is 


Eapfpp  - m £x  - - 0x3  - -<P)£)p  * £t  ' 


(7.10) 


where  (see  §3) 


r°°  i 

P = J E^y  (s)  y (0 ) ) ds . 


(7.11) 


The  equilibrium  density  is  (8  = kT) 


f = exp 
eq 


[-»(£  + + nr)" 


(7.12) 


We  require  that  f be  a solution  of  (7.10)  and  obtain  (the 

eq 

"fluctuation-dissipation"  result) : 


e ap  = 2kTy. 


(7.12) 


The  steady  states  of  the  averaged  equations  (7,8,9)  are 


where  we  have  made  the  assumption  that  k > 0,  a > 0. 

We  assume  that  when  n = nc»  k(nc)  = 0.  Then  (0,0)  is  a critical 
type  steady  state.  At  n = nc,  we  have  a critical  harmonic 
oscillator. 

kT  i 

We  now  nondimensionalize  (7.9).  We  let  e = — <<1  be  a 

En 


small  parameter,  where  EQ  is  some  reference  energy.  Introducing 
dimensionless  variables  by 


In  the  sequel,  we  drop  the  primes  in  (7.16).  Since  k'  ® k,  when 


n = nc,  k ' (nc)  = 0. 

We  now  seek  a solution  of  (7.16)  of  the  form 

,2 


f(t,x,v)  = exp 


^(S»4-=S-- 


( x , t ) (7.17) 


where  4*, a,  6 and  zn  are  to  be  determined.  Following  the  procedure 
in  §5,  we  obtain 

4*  + vf  - (k(r|)x  + ax'*  + yv)  41 


+ y(V3  - a'V  - B)4»J  = 0, 


(7.18) 


Let  us  now  specialize  to  n = n , k = 0;  i.e.  the  critical 

c 

harmonic  oscillator.  Then  a = 6 = 0 in  (7. 17) and  (7.18).  The 
ray  equations  become 


dx 

dt 


= v 


§£  = -ax3  - yv  + 2 y 4'v 4* 3 


Slit  = v4<  + f + \j» 

dt  x dt  v Yt 


3F1  * -6  Vv'’2  ' 3&t\ 


d'V 

ar-  = -"v*2  - Tx  + *Tv 


(7.19) 


L. 
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By  integrating  the  ray  equations  from  an  initial  point  x(tQ)  = xQ , 
v(tQ)  = vQ/  we  obtain  the  conditional  density  f (x,  t,v  ,Xq  ,v  • ) . Then 
following  the  procedure  in  §6,  we  can  obtain  the  correlation 
function.  Our  solution  thus  allows  the  calculation  of  correlation 
functions  at  critical  points.  Once  the  correlation  function  is 
known,  we  obtain  the  spectrum  of  the  oscillator  by  Fourier 
transform. 
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APPENDIX  A 


MARGINAL  AND  CRITICAL  TYPE  DYNAMICAL  SYSTEMS 

In  this  appendix,  we  give  exact  conditions  for  marginal  and 
critical  type  dynamical  systems.  Our  work  generalizes  the  scheme 
of  Kubo  et.  al. (1973) . 

MARGINAL  TYPE  DYNAMICAL  SYSTEMS 

The  deterministic  evolution  of  the  macrovariables  is 
governed  by 

x = b (x,  n)  (A.  1) 

where  n R1  is  a parameter.  Equation  (A.l)  may  have  three  steady 
states,  QQ  (n)  , Q^n)  and  P2>  Let  Bk  be  the  matrix  (b1,^) 
evaluated  at  QQ , Q1  or  P2  (k  = 0,1,2).  We  assume  that: 

• For  all  values  of  n,  B2  has  two  real  negative  eigenvalues. 
Although  P2  may  depend  upon  n,  P2  is  always  bounded  away  from  the 
other  steady  states. 

• As  q I 0 , the  distance  between  Qq  (h)  and  Q^(n)  decreases. 

When  n = 0,  Qq  and  coalesce  and  annihilate  each  other  (i.e. 

when  n < 0 , (A.l)  has  one  real  and  two  complex  steady  states). 

• When  n >0,  has  two  real  negative  eigenvalues  and 
has  one  real  positive  and  one  real  negative  eigenvalue.  When 

n = 0,  bp  = B^  has  one  zero  and  one  real  negative  eigenvalue.  The 


r 


eigenvector  corresponding  to  the  negative  eigenvalue  has  positive 
slope.  The  double  point  Qq(0)/Qj,(0)  is  called  a saddle  node  (3). 

A deterministic  system  satisfying  the  above  assumptions  will 

be  structurally  similar  to  the  system  sketched  in  Fig.  2. 

The  above  conditions  can  be  reformulated  by  a change  of 

coordinates.  Define  the  y^  axis  in  the  direction  of  the  eigen- 

2 

vector  of  the  non-negative  eigenvalue  of  . The  y axis  is  in 
the  direction  of  the  eigenvector  of  the  negative  eigenvalue  of 
B^,  with  the  origin  at  Q^.  Then 

y = b (y, n)  (A. 2) 

is  the  deterministic  system  in  the  new  coordinates.  The  system  is 
of  the  marginal  type  if: 

1)  det  (b1 , j (Qr0)  ) = 0 

2)  b 1,1(Q1,0)  =b  2,1(Q1,0)  = 0 (A. 3) 

3)  b 2 * 2 (Q]_ ' 0 ) * 0 

4)  b 1,11(Q1,0)  - b 2,11(Q1,0)  ? 0. 

The  conditions  (A. 3)  have  the  following  interpretation.  Condition 
1)  indicates  that  the  original  system  has  a zero  eigenvalue. 
Condition  2)  indicates  that  when  n = 0 the  linear  dynamics  in  the 


A-2 


y direction  vanish,  condition  4)  indicates  that  these  dynamics 
are  quadratic.  Condition  3)  indicates  that  the  second  eigenvalue 
is  non-zero. 


CRITICAL  TYPE  DYNAMICAL  SYSTEMS 


The  macrovariables  evolve  according  to  a deterministic 


kinetic  equation 


x = b (x,  n#  <5) 


(A. 4) 


where  n#  6 are  one  dimensional  parameters.  The  entire  bifurcation 
set  of  equation  (A. 4)  is  still  unknown  (3).  The 
physical  systems  of  interest  here  motivate  the  following 
assumptions : 

• For  some  values  of  n > S,  (A. 4)  has  three  steady  states 
PQ(rif5),  P^fn/S)  and  P2(n,<S). 

If  B^  = (b1^  ) evaluated  at  P^,  then  when 
the  three  steady  states  are  distinct,  BQ  and  B2  have  real  negative 
eigenvalues.  B^  has  one  real  negative  and  one  real  positive 
eigenvalue.  The  eigenvector  corresponding  to  the  negative  eigen- 
value has  positive  slope. 

• As  r)/  6 vary,  two  of  the  steady  states  may  coalesce  and 
annihilate  each  other.  This  behavior  is  analogous  to  the 
marginal  bifurcation. 


A- 3 


I 

1 


• As  Hi  <5  vary,  all  three  steady  states  may  move  together 
and  coalesce  when  n = <5  = 0.  At  the  critical  bifurcation, 

= (b^J  has  a zero  eigenvalue.  We  assume  that  the  steady 
state  remaining  after  the  critical  bifurcation  is  a stable  steady 
state. 

A deterministic  system  satisfying  the  above  postulates  will 

be  structurally  similar  to  the  one  sketched  in  Fig.  3. 

The  above  properties  can  be  restated  in  terms  of  a new 

coordinate  system  as  follows.  The  y'*'  axis  is  in  the  direction  of 

2 

the  eigenvector  of  the  non  negative  eigenvalue  of  . The  y 
axis  is  in  the  direction  of  the  eigenvector  of  the  negative 
eigenvalue,  with  the  origin  at  P^.  The  deterministic  evolution 
is  then 


y = b (y , nr  S)  • 

A dynamical  system  is  a critical  type  system  if: 


1)  det(b  ,j(P1,0,0))  = 0 


2)  b 1,1(P1,0,0)  = & 2 3 4,1(P1,0,0) 


= b 1,11(P1,0,0)  = b 2,11(P1,0,0) 


3)  b 2,2(P1,0,0)  ft  0 

1 y 2 

4)  b '111  “ b '111  * °* 


(A. 5) 


(A. 6) 


These  conditions  have  the  following  interpretation  1)  indicates 
that  the  system  has  a zero  eigenvalue,  while  condition  3)  indicates 


that  the  second  eigenvalue  is  non-zero.  Condition  2)  indicates 

1 


that  the  linear  and  quadratic  dynamics  in  the  y direction  vanish, 
while  4)  indicates  that  the  dynamics  are  cubic. 
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